Phase transition in the collisionless damping regime for wave-particle interaction 
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Gibbs statistical mechanics is derived for the Hamiltonian 
system coupling self-consistently a wave to N particles. This 
identifies Landau damping with a regime where a second or- 
der phase transition occurs. For nonequilibrium initial data 
with warm particles, a critical initial wave intensity is found : 
above it, thermodynamics predicts a finite wave amplitude in 
the limit N — > oo ; below it, the equilibrium amplitude van- 
ishes. Simulations support these predictions providing new 
insight on the long-time nonlinear fate of the wave due to 
Landau damping in plasmas. 
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Landau damping is a striking phenomenon first evi- 
denced in plasma physics, but also lately in various sys- 
tems such as nonlinearly coupled oscillators |l]] or neu- 
tral rarefied gases (2|. In plasma physics, its classical 
treatment involves a continuous media analysis through 
Vlasov-Poisson equations. It was shown by O'Neil that 
(near-)resonant particles play a major role in the time- 
asymptotic, generically nonlinear regime experienced by 
kinetic systems submitted to Landau damping ||. This 
has been debated recently by several works. In partic- 
ular, Isichenko argued for algebraic final decay of the 
electric field in the nonlinear regime 0, under the cru- 
cial assumption of an asymptotically vanishing field 
though in numerical investigations by Manfredi || the 
electric field was seen to damp exponentially in the linear 
regime and then to oscillate around a finite amplitude in- 
dicating a Bernstein-Greene-Kruskal (BGK) equilibrium. 
Lately, Lancellotti and Doming |7]] have shown that both 
types of time-asymptotic states were eligible depending 
on initial conditions. 

In this Letter we address this issue from a differ- 
ent viewpoint and relate both possible time-asymptotic 
states (with zero or finite field) to different phases in a 
phase transition picture. Rather than using the large 
number of weakly interacting particles N to express the 
dynamics by a kinetic equation for a distribution function 
we note that the intrinsic chaos in this dynam- 
ics ensures an effective ergodicity allowing a Gibbsian 
equilibrium statistical mechanics treatment. This analy- 
sis provides further insight on the results of the kinetic 
approach. 

As real plasmas are made of a finite number of par- 
ticles and as these graininess properties, and notably 



the existence of spontaneous emission, are discarded by 
their classical vlasovian treatment, a model taking these 
properties into account has been derived ||,[|- The first 
steps in this direction were made by O'Neil et al. pQ ] 
in order to treat the cold beam-plasma instability, and 
a Hamiltonian formalism was first introduced by Myn- 
ick and Kaufman pd[ . This Hamiltonian system mod- 
els the interaction between M (bulk) Langmuir waves 
and N quasiresonant (tail) particles and has already pro- 
vided insight into linear Landau growth and damping, 
the cold beam-plasma instability |12| and Van Kampcn 
modes |J. From the physical viewpoint, this approach 
has the advantage to treat the waves as dynamical objects 
rather than subordinate phenomena, i.e. excitations of 
the bulk plasma. In the kinetic limit N — > oo, the coupled 
(many-body) wave-particle dynamics approaches the cor- 
responding dynamics generated by the Vlasov equation 
over any finite time interval ]l3| . 

For simplicity, we concentrate on the single-mode case 
(M = 1), which is paradigmatic to study Landau damp- 
ing. Besides, this case is of special interest for devices 
such as traveling wave tubes where a single mode can 
be selected For N identical (tail) particles moving 
on the interval of length L = 2n with periodic boundary 
conditions T := R/(2vrZ), with unit mass and charge, 
and respectively position x r and momentum p ri inter- 
acting with one wave of natural frequency luq, unit wave 
number, phase 9 and intensity /, the Hamiltonian reads 
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The first term describes ballistic motion of the particles, 
the second the oscillations of the free wave (harmonic 
oscillator) and the third couples them; here r\ is a small 
parameter denoting the ratio of the tail density over the 
plasma density. The bulk of the plasma is considered 
as a linear dielectric supporting plasma waves. The wave 
frequency loq (which is the plasma frequency ui p ) provides 
the natural reference time scale, which can be fixed with 
no loss of generality. Thus ujq 1 will be the time unit, so 
that the phase velocity of the free wave is also 1 . 

The self-consistent system can be viewed as a set of 
nonlinear oscillators coupled through a mean-field, as the 
equation of motion for any particle r obtained from (|l|) 
reads 



= -y/2rjI/N sin(a 
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i.e. the equation of a pendulum in the field with strength 
y/2r)I/N and angle 0. The wave evolves consistently 
with the particles. For this mean-field model, the ki- 
netic vlasovian limit ]T^ ] corresponds to taking the limit 
N — > oo while the acceleration felt by any particle should 
remain finite. This induces the natural kinetic scaling for 
the wave intensity I = O(N) and gives its kinetic equiv- 
alent as the intensive quantity ip := I/N. 

Most importantly, this approach enables an equilib- 
rium statistical mechanics treatment, provided the dy- 
namics is effectively ergodic. In this Letter we derive the 
equilibrium Gibbs predictions for one field's observable, 
its intensity I. This study predicts a second-order phase 
transition in the regime corresponding to damping. We 
exhibit a critical initial field intensity, above which par- 
ticle trapping induces a finite field amplitude, whereas 
below it, in the kinetic limit of an infinite number of par- 
ticles, the asymptotic field would vanish. These results 
complement the dynamical ones [§,0, and our estima- 
tions are fully supported by numerical simulations based 
on a fourth-order symplectic scheme Q . 

The Liouville sub-manifold of the phase space is fixed 
by the two usual constants of the motion, total energy 
E = H and total momentum P = 53r=iPr + I- To 
simplify calculations, we turn to the canonical ensemble, 
relying on the equivalence of ensembles when systems get 
large |l6[ ], and introduce the temperature T (ks = 1) of 
the tail particles-wave system. The canonical measure 
(constrained by constant momentum P) reads 



d/j, c 



-H/T 5 |p 



N 



N 



2_.Pr ~ I I dldO J^J dp r dx r (3) 



where the phase space variables (p, /, x, 9) evolve in A 



xFxT. With the intensive variable a := P/N, 



the canonical partition function Z C (T, a, N) reads 
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where /(ct,T» := {a-^f /2 + uj^j -T\nl a {^2^ /T) ; 
I n denotes the modified Bessel function of order n. For 
large TV, with tp :— y/2rjip/T, one finds the minimum 
of g(<r,T,tp) := f[a, T,Ttp 2 /(2rj)}, by solving the saddle- 
point equation 



d v g = tpT 
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The minimum of g is reached at tp — iff a < luq and 
T > T r , where 



V 



2 (u) - a) ' 



(6) 



Otherwise the minimum is reached for a certain 
(p*(a,T) > 0, defined implicitly by equation (||). The 
canonical averages for the wave intensity are respectively 
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Consequently, when a < ujq, a second order phase transi- 
tion p6| , p~7j occurs at T c , with order parameter ip = I/N, 
i.e. the suitably normalized intensity of the wave [ p^| . 
Above T c , the intensity of the wave loses its extensivity : 
the particles and the wave decouple in the sense that the 
wave controls just one among N + 1 degrees of freedom. 
On the contrary, for T < T c , the wave is macroscopic 
and large enough to trap particles. When T — > T c , (0) 
exhibits a formal divergence signaling the crossover from 
intensive to extensive scaling (|§|). 

Finally we express the equation of state of the 
particles-wave system in terms of the Gibbs average of 



the energy density h 
Thus, for tp = 0, 



H/N , using (H) c = T 2 9 T lnZ c 



WcL = o = ^ 2 /2 + T/2, 



(9) 



where cr 2 /2 is interpreted as the mean centre-of-mass ki- 
netic energy of the particles- wave system and T/2 as the 
thermal agitation energy in this frame. For tp* ^ 0, in 
the limit of large TV, one obtains similarly 



(h) 



g[tp*(a,T)]+T/2. 
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Before considering the implications of our results, let 
us discuss the conditions under which they hold. The 
Gibbs Ansatz assumes that particles explore the avail- 
able (x,p) space thanks to the nonintegrability of their 
dynamics. For physical applications, it suffices that they 
explore a 'large' part of the (x,p) space ; in our case, 
that they be able to wander, say, within a velocity range 
of two standard deviations on either side of their aver- 
age velocity. The constraint J2 r Pr + I — P ensures that 
the kinetic energy contribution to (|^) is a Maxwell distri- 
bution centered on a — tp, with standard deviation T 1//2 
for p r . One therefore expects some good mixing prop- 
erties, due to a trapping/detrapping mechanism [ p2| , if 
T 1 / 2 is of the order of the velocity width of resonance 
(^), namely (r/ip) 1 ^. This equivalently corresponds to tp 
being of order 1. In discussing numerical simulations, we 
shall note that relaxation towards equilibrium occurs for 
initial data roughly realizing this condition. 

Now we turn to the manifestation of this phase transi- 
tion in the wave-particle dynamics. The problem raised 
by Landau damping is an initial value problem (hence 
a nonequilibrium one), where a packet of particles, spa- 
tially homogeneous with velocities distributed according 
to a certain fo(v), interacts with a wave launched with 
intensity ipo > 0. This determines E/N — (v 2 )o/2+u)o4>o 
and <7 = (v)q + ipQ. The condition a < lu () implies 
(v) < uj , which corresponds typically to a bunch of 
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particles, centered on the phase velocity of the wave, ini- 
tially distributed according to a decreasing fo(v) and thus 
inducing Landau damping on the wave. The initial par- 
ticles temperature is T := (v 2 )o — (v)q. 

Given E/N and a, the equation of state (||) in the 'high 
temperature' regime expresses the equilibrium tempera- 
ture T in terms of initial data 

T = T + 2^o(^o ~ (v)o) - (11) 

Using (^), defining a := Wo — {v)o > w o — & > and 
putting a± :— a±y/a 2 + To, condition T > T c reduces to 
2(ipo — a)(ipo — a+)(ipQ — a^) > 77. Therefore, provided the 
(near-resonant, tail) particle distribution is warm enough 
(2aTo > 77), if ipo < "00c in the limit of large TV, the 
equilibrium field (ip) c vanishes, whereas for ipo > V'Oc the 
equilibrium field remains finite. The critical initial wave 
amplitude "00c is defined by 

2 (ip 0c - a) (tpoc - a+) (ip 0c - a_) = 77. (12) 

Moreover, whatever the initial conditions, one deduces 
that ii ipo > a then the asymptotic field will always be 
finite. 

These results are illustrated for an initial normalized, 
warm particle distribution fo(v) = C\ — Cij{2 — v), 
with C\,C2 > 0, over the range Vi < v < V2 (with 
t>2 — loq = ^0 — v i) and fa(v) — outside of it. The 
linear Landau rate || for amplitude \ftp~ — \/I/N is 
7l = (tt /2)LUor)fo(uJo)- In our dynamical simulations, 
w = 1, 17 = 5.02 • 1(T 4 , Vi = 0.75, v 2 = 1.25 and 
7l = —1/200 (from which C\ and Ci are deduced). Then 
(v) = .931 and (v 2 ) = .882, thus T = 1.50 • 10~ 2 . The 
resulting critical mean intensity is t^oc = 5.59 ■ 10~ 2 . 

Running simulations for different values of ipo and N 
yields the plot of Fig. [l] for asymptotic values of ip. 
This figure shows a clear agreement with thermodynam- 
ics predictions in the kinetic limit N — > 00 for both 
quasi-ballistic and trapping regimes. In the critical re- 
gion where ipo ~ V'Qc (dashed lines), one observes strong 
metastability effects and critical slowing-down, whereas 
additional finite- N corrections smooth the transition : in 
this region, initial field intensities are far larger than their 
equilibrium averages. Moreover, for the range of ipo con- 
sidered on Fig. 0, the analytic expression for T happens 
to differ only slightly from To. Therefore the 'mixing' 
condition corresponds to ^/2rjipo/To being of order 1. 

This mixing condition is roughly fulfilled for the runs 
of Fig. ^| displaying the temporal evolution of In \fip and 
progress towards equipartition for different values of N . 
Along these runs, the relative energy variation is less than 
2 • 10" 7 . On the inset, after the initial linear Landau 
damping at the prescribed rate, one observes trapping 
oscillations that could suggest the setting up of a BGK 
equilibrium. However, on longer times, the sweeping of 
the phase space is sufficient to drive the system to a far 
different, 'high' temperature, regime where I is no longer 



extensive. Critical slowing-down reflects on the relax- 
ation times diverging with N . Analogous pathological 
relaxation properties in the limit N -— > 00 have been re- 
ported for a similar mean-field model |lj] in the subcrit- 
ical energy domain where chaos is paradoxically shown 
to be maximal [ |l8| , |l9| ]. 

It is surprising to note the good ergodic behavior of ip 
for the initial domain where ipo <C ipo c , where no mixing 
can be expected as the system is close to integrability. 
Here one recovers the effective linear Landau damping 
and ergodicity is ensured dynamically by a perturbative 
analysis P,|20|]. 

In the trapping regime, the equilibrium level is quickly 
reached as a substantial fraction of particles is liable to 
undergo separatrix crossings J|2| . Nevertheless ipo is lim- 
ited by the velocity range chosen for the tail particles 
p3| . One checks that the trajectories corresponding to 
the velocity borders of our distribution in the Boltzmann 
(x,p) space correspond to KAM tori for the effective 1.5 
degrees of freedom Hamiltonian dynamics (J2J) with time 
dependent parameters {0,1) |2l|| . 

Finally, denoting by tbo = (^rjipo) -1 ^ 4 the trapping 
time corresponding to the initial mean wave intensity and 
putting q — |7l|tbo, it turns out that the phase transi- 
tion occurs here for q — q c ~ 0.06. This agrees with 
the qualitative dynamical threshold drawn (in a non self- 
consistent way) by O'Neil ||] separating a regime domi- 
nated by trapping (if q < q c ) from a regime where linear 
predictions are effective (if q > q c ). Existing numeri- 
cal self-consistent simulations of the damping of a single 
wave are limited to the early stage of wave evolution and 
would typically infer a larger q c (e.g. q c — 0.77 for a 
linear /o(i>) in |24|]). As clearly accounted for by Fig. |2| 
and inset (for which q = 0.14) , our new approach relates 
this apparent discrepancy to critical slowing-down. Our 
intuition is that an algebraic decay of the order param- 
eter ip (in the limits — > 00 and t — > 00) may be a 
dynamical signature of the vicinity of critical point, as in 
Landau second order phase transitions' picture. In a pa- 
rameter regime q i^0.06, Isichenko's conclusions Q may 
then apply here. 

In conclusion, using the self-consistent wave-particle 
Hamiltonian (|l|) to describe this system has revealed a 
new phenomenon, the phase transition associated with 
the Landau damping regime. This phenomenon eludes 
the usual Vlasov-Poisson description, which actually 
treats all particles on the same footing, as a Coulomb 
system. But the wave-particle interaction is effectively 
very different for near-resonant particles (trapping time 
scale (277-0) -1 / 4 and Landau linear scale |7l| _1 ) and for 
bulk particles (adiabatic or non-resonant averaging time 
scales, exponentially longer than resonant ones): our self- 
consistent Hamiltonian describes the system on the phys- 
ical time scales for the wave evolution M . 

Moreover, our finite-A approach shows that, in a real 
plasma, the wave will never damp completely but will 
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eventually fluctuate around a finite TV-dependent thermal 
level due to spontaneous wave emission. 
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the resonant region to extend as far as the bulk of the 
plasma, as this would also break the validity of the self- 
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FIG. 1. Time averages of normalized intensity ip reached 
for long times versus i/'o for TV = 16000 (triangles), TV = 32000 
(squares) and TV = 64000 (stars). Filled symbols for simula- 
tions times of 2 ■ lO 8 ^" 1 . Open symbols for run times of 
5 • 10 Wo . Canonical estimations are drawn in continuous 
line for TV = 16000, 32000 and 64000 particles. Near the 
transition, dashes indicate that estimate (Q) is not expected 
to be accurate. 



FIG. 2. Time evolution of normalized amplitude 

v 7 ^ = \/T/N for TV = 16000 (dots), TV = 32000 (thick line) 
and TV = 64000 (thin line), for Vo = 1.56- 10" 3 < ip 0c . Inset : 
initial evolution, including brief initial linear damping stage. 
On longer times, for TV = 32000, ip wanders around half the 
value associated to TV = 16000, in agreement with canonical 
predictions. Relaxation towards equilibrium is much longer 
for TV = 64000 due to critical slowing-down near the phase 
transition. 



5 



G 



